﻿using LatoolNet;

namespace NumericalSolution {
    class Diffusion1DCN : Diffusion1D {
        public override double[] Next() {
            base.Age();

            int rownum = p_numberOfNodes;
            int colnum = p_numberOfNodes;
            int bandWidth = 3;
            Matrix mat = new Matrix(rownum, colnum, bandWidth);
            Matrix vec = new Matrix(colnum, 1);

            double s = (p_alpha * p_deltaTime) / (p_deltaX * p_deltaX);

            mat[0, 0] = 1;
            vec[0, 0] = 0.0;
            for (int i = 1; i < p_numberOfNodes - 1; i++) {
                mat[i, i - 1] = -s;
                mat[i, i] = (2 + 2 * s);
                mat[i, i + 1] = -s;
                vec[i, 0] = s * p_previous_u[i - 1] +
                            (2 - 2 * s) * p_previous_u[i] +
                            s * p_previous_u[i + 1];
            }
            mat[p_numberOfNodes - 1, p_numberOfNodes - 1] = 1;
            vec[p_numberOfNodes - 1, 0] = 0.0;

            //LUFactorization.Solve(mat, vec);
            LinearEquation.Solve(mat, vec);

            for (int i = 0; i < p_numberOfNodes; i++) {
                p_current_u[i] = vec[i, 0];
            }

            mat.Dispose();
            vec.Dispose();
            return p_current_u;

        }
    }
}
